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Realistic black hole collisions result in a rapidly rotating 
Kerr hole, but simulations to date have focused on nonro- 
tating final holes. Using a new solution of the Einstein ini- 
tial value equations we present here waveforms and radiation 
for an axisymmetric Kerr-hole-forming collision starting from 
small initial separation (the "close limit" approximation) of 
two identical rotating holes. Several new features are present 
in the results: (i) In the limit of small separation, the wave- 
form is linear (not quadratic) in the separation, (ii) The wave- 
forms show damped oscillations mixing quasinormal ringing 
of different multipoles. 
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The study of black hole collisions with numerical rel- 
ativity will give unprecedented insights into the non- 
linear workings of relativistic gravitation in strong field 
situations and is expected to provide important infor- 
mation about the waves that may be observable with the 
gravitational wave detectors currently under construction 
p|. In the difficult task of developing collision codes it 
has been found useful to compare numerical relativity re- 
sults with those of "close limit" calculations || . In these 
calculations the black holes are taken to start at a small 
separation, and the initial separation is used as a per- 
turbation parameter. The collisions of primary interest 
are those with high angular momentum that end in the 
formation of a rapidly rotating Kerr black hole. This is 
based on the expectation that such collisions are the end- 
point of binary inspiral and that collisions with high an- 
gular momentum will produce much more gravitational 
radiation than collisions leading to a nonrotating hole. 
But the close limit method has so far been applied only 
to collisions resulting in a nonrotating, or slowly rotating 
hole. There are two reasons for this. 

The first reason is a technical one. The starting point 
of the dynamics of a collapse is a spatial metric ( 3 ) gij and 
an extrinsic curvature K+j that solve Einstein's initial 
value equations. For the close limit perturbation method 
to be used for rotating holes, these initial value solutions 
must be a perturbation of a Kerr spacetime. But the 
prevalent method for specifying initial data, the confor- 
mally flat prescription of Bowen and York Q , is incom- 
patible with the Kerr solution, for which the standard 
slicing, at least, is not conformally flat. A family of close 
limit perturbations of Kerr initial data, then, could not 



be given in the Bowen- York formalism. 

This technical difficulty is no longer a barrier. New 
axisymmetric initial data solutions representing two ro- 
tating holes have recently been given by Baker and Puzio 
H , and by us |g] . In both types of solutions close limit se- 
quences of initial data can be constructed for holes start- 
ing a head-on axisymmetric collision. More specifically, 
a sequence can be constructed for different initial separa- 
tions, with the close limit of the sequence being the Kerr 
spatial geometry and extrinsic curvature outside a Kerr 
horizon. We have chosen to start with the close limit ini- 
tial data set of Ref. || based on the Brill-Brandt-Seidel 
form of the 3-metric and of the nonvanishing compo- 
nents of extrinsic curvature, 
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in which q, Hp, Hp are functions of r, 9. For the ap- 
propriate choices ^KjQk,Hek,Hfk, these become the 
metric and extrinsic curvature for a slice of Kerr space- 
time at constant Boyer-Lindquist || time. In this case 
the 9 coordinate of (|l|) and (Q) is the same as the Boyer- 
Lindquist 9, and the Boyer-Lindquist radial coordinate 
r is a rescaling of f, depending on the Kerr parameters 
M, a, as given in Ref. |J. 

In our approach to an initial value solution we take 
q, He, Hp to have precisely their Kerr forms qn, 
Hek,Hfk- ft turns out that this guarantees that all 
of the initial value equations are satisfied except one, 
the Hamiltonian constraint, an elliptic equation for $. 
Boundary conditions for this equation must be chosen in 
order to complete the specification of a solution. One 
condition is asymptotic flatness. If for the other con- 
dition we choose to have <E> ~ K S j ng /f as f — > 0, then 
we find the Kerr solution if K S i ng is chosen appro- 
priately. (Otherwise we find a distorted Kerr hole.) To 
arrive at an initial value solution representing two initial 
holes, we choose two point-like singularities on the z axis 
in the f, 9, plane, and we place them at z = ±zq (that 
is, at f = Zo,6 = 0,7r). The parameter Zq, then, de- 
scribes the initial separation between the two holes. It is 
shown in Ref. |J that the zq — > limit of the initial solu- 
tion is the Kerr geometry, outside the Kerr horizon, with 
mass M and spin parameter a (the parameters chosen 
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for the functions qk,Hek,Hfk)- In the limit of small 
zq, the 3-geometry and extrinsic curvatures are pertur- 
bations away from the Kerr forms. This provides the 
starting point for evolving perturbations on a Kerr back- 
ground. 

The second reason that close limit calculations 
have previously been limited to collisions forming a 
Schwarzschild (nonrotating) hole is that no numerical 
relativity results were available of collisions leading to 
Kerr formation. This can, of course, be viewed as an 
opportunity to provide the only available - though very 
limited - answers about the nature of collisions to form 
Kerr holes. But our real motivation in presenting these 
results is to continue the useful interaction of perturba- 
tion studies and numerical relativity. In particular, we 
hope that these results stimulate work on fully nonlinear 
evolution of initial data solutions amenable to close limit 
perturbation methods, and hence to comparisons. 

Linearized evolution for perturbations of the Kerr ge- 
ometry is carried out with the Tcukolsky equation 
rather than the Zerilli or Regge- Wheeler equa- 
tion of Schwarzschild perturbations. The introduction 
of the Teukolsky equation entails two new elements in 
a calculation. Since the Kerr background is not spheri- 
cally symmetric, its perturbations cannot be decomposed 
into spherical harmonics Jl2]] . The dependence on a po- 
lar angle 9 cannot therefore be separated. The linearized 
numerical evolution must be carried out as a solution of 
a 2+1 variable (r, 9, t) linear hyperbolic differential equa- 
tion, rather than a 1+1 variable (r, t) equation for each 
multipolc of spherical perturbations. (Since we are dis- 
cussing the evolution of perturbations of the Kerr space- 
time, here and below the coordinates r, t, 9 are the stan- 
dard, i.e., Boyer-Lindquist || coordinates for the Kerr 
spacetime.) To solve the 2+1 partial differential equation 
we use an existing and thoroughly verified 2+1 evolution 
code fl3fl . Though the computational requirements of 
this code are much greater than for 1+1 codes, they are 
handled easily by the RAM and speed of large modern 
workstations. It should be noted that the essential diffi- 
culties of numerical relativity lie in the nonlinear equa- 
tions of Einstein's theory. None of the worst difficulties 
affect our linearized version of the evolution equations. 

The second new feature to be faced is the way in 
which initial conditions must be handled. In the Zerilli 
|Tof or Regge- Wheeler [[ll) formalism, evolution is car- 
ried out for a "wavefunction" that is defined in terms 
of metric perturbations. It is relatively straightforward 
to compute the initial wavefunction and its initial time 
derivative from an initial value solution for ( 3 ' gij , . 
In the Teukolsky formalism the quantity that is evolved 
with a wave equation has a somewhat different charac- 
ter. This quantity originates in the Newman-Penrose |Ti| ] 
formalism that encodes information about curvature into 
five complex fields ipQ, ... ,1/14 that are projections of the 
Weyl tensor (equivalent to the Riemann curvature tensor 



Ra.p-y5 in vacuum) onto a null tetrad of complex vectors 
^,71^, , fh^. Of particular importance here, is the pro- 
jection 

V>4 = R a p 1 s'n a fh^n 1 fh s , (3) 

since the Teukolsky formalism evolves the wave function 

* = (r-ia cos 60 4 i/>4 ■ (4) 

It turns out that ipi (and hence ^) is invariant to first 
order under first order perturbations of the tetrad, so the 
specification of initial data amounts to the computation 
of the initial value of R a f3~fS, and its initial time deriva- 
tive. The computation of the initial "J is carried out using 
the Gaufi-Codazzi equations Jl5| which relate the Rie- 
mann curvature to ^ gij, Kij of the initial value solution. 
Finding the initial requires that we know the initial 
time derivative of the Weyl tensor. To find this we substi- 
tute the first order vacuum Einstein evolution equations 
for the 3-metric and the extrinsic curvature into the ex- 
plicit formula for dt% ■ Though straightforward in prin- 
ciple, the calculation of the initial \f r and of dt^f from the 
initial ^ gij, Kij is very lengthy, and subject to error. To 
carry out the computation, a Maple script was written 
to find W from the initial data functions &,q, He, Hp. 
The output of the script was checked against analytic 
expressions known for a = ]l6fl . For a ^ the result 
was checked numerically against a very different Maple 
script, kindly supplied by Gaurav Khanna, in which ip^ 
is computed from the linearized perturbations on a Kerr 
background jl7| . In both comparisons the agreement was 
well within the expected numerical accuracy. 

For values restricted to the region outside the event 
horizon of the Kerr background, it is useful to intro- 
duce the radial variable r* defined Jl8| in terms of the 
standard (Boyer-Lindquist) radial variable r by dr* jdr = 
(r 2 + a 2 )/(r 2 - 2Mr + a 2 ). The r* variable agrees with r 
in the limit of large r and approaches —00 at the unper- 
turbed horizon r + = M + v 'M 2 — a 2 . Since the value of 
^ diverges as r 3 for large r we introduce Q = Mvf/r 3 . 

We have chosen to evolve initial solutions represent- 
ing holes of equal mass and spin, so that each collision 
is characterized by the set of parameters a, M, zq (where 
M and a are the mass and spin parameters of the fi- 
nal Kerr hole), or by the two dimensionless parameters 
a/M,zo/M. Since ^ and dt^f are complex, there are 
four functions of r, 9 contained in them. Due to the sym- 
metry of the collision these functions all have symmetry 
properties with respect to the 9 = tt/2 equatorial plane. 
Specifically, Re\E , ,Re9t^' are symmetric with respect to 
the equatorial plane, while Im^Imc^ are antisymmet- 
ric. The plot of Im<9 f <2 as a function of r* , in Fig.[l], 
shows the antisymmetry of Im^^ and the vanishing of 
'J at the horizon and at large r. 

In the case a = the conformal factor <E> can be written 
in analytic form, and it can be shown that in the zq — > 
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limit, the relative deviation from the Schwarzschild ge- 
ometry is quadratic in zq. In the numerically computed 
solutions for $, we find that the relative deviation 5 is 
generally linear in zq. Since this dependence must revert 
to a quadratic one for a — > 0, we infer 



S = F x {r, B- a)z Q + F 2 (r, 9; a)z 2 + 0{4) 



(5) 



where F\ vanishes, but F 2 does not, as a — > 0. A numer- 
ical confirmation of this relationship is shown in Fig. ^. 

The values of \& and as described above are next 
used as Cauchy data for the Teukolsky equation, and the 
equation is solved numerically to find Q(r, 8, t). The real 
and imaginary parts of \& correspond to the two distinct 
linear polarization modes of gravitational waves, and 
both ReQ (solid curve) and ImQ (solid curve) are shown 
in Fig.||. A semilog plot is used to show most clearly 
the nature of the oscillations and exponential damping 
of the waves, after an initial transient. Some features 
of Fig. H have a simple explanation: Although the Kerr 
background is not spherically symmetric, perturbations 
at a single frequency can be decomposed into angular 
functions and one can find the complex quasinormal fre- 
quencies for each angular function. For ReQ both the 
damping (slope of the line connecting the peaks) and the 
oscillation rate (spacing of zeroes) agree to within around 
1% with the frequency of the least damped quasinormal 
frequency @Wj = (0.388 + i0.086)/M for the £ = 2 ax- 
isymmetric mode of the a/M = 0.6 hole that is formed in 
the collision. The explanation of ImQ is less simple, but 
more interesting. It has the appearance of a superposi- 
tion of ringing at oj 2 and at ui 3 = (0.618 + z0.088)/M. For 
t/M > 50 the dashed curve in fact is indistinguishable 
from the function 

4.39 x 10" 4 sin(0.388t/Af + 0.463) exp(-0.086t/Af) 
+ 8.26 x 10~ 4 sin(0.618t/Af - 1.951) exp(-0.088i/M) . 

This superposition comes about because the Teukolsky 
equation mixes real and imaginary parts, through a term 
involving i acosOdt^ . The angular dependence of the 
imaginary part of the initial data gives rise to small ImQ 
evolution characterized by W3, but during evolution, the 
oscillations of ReQ, at LO2 are fed into ImQ by the Teukol- 
sky equation. The reverse mixing of ImQ into ReQ is not 
of importance because ImQ is much smaller than ReQ, 
and mixing is weak. 

Figure || gives the computed total gravitational wave 
energy contained in the waveforms, as a function of zq 
for spin parameter values a/M = 0,0.1,0.2,0.4. The 
a = problem is equivalent to computations using 
Schwarzschild perturbation methods, and the results of 
those methods |2(J agree with the a = points in Fig.|] 
to better than 1.5%. The figure is somewhat misleading 
regarding the dependence of radiated energy on angu- 
lar momentum. The variable zq is a formal measure of 



separation, not directly connected with a physical sep- 
aration. A physical separation, for example the proper 
distance between the intersection of the axis and the ap- 
parent horizon, depends on both zq and a. Figure ||, 
therefore, does not show the influence of rotation on ra- 
diation, but it does give results for well specified initial 
data, results that could be evolved with numerical rel- 
ativity on supercomputers. We hope that comparisons 
with such results will be available soon. Another appli- 
cation of these results will be a comparison of predicted 
waveforms and energies using the slow-rotation approx- 
imation (i.e., treating rotation as well as separation as 
a perturbation), a comparison that will reveal whether 
slow-rotation methods can be used for fairly rapid ro- 
tation. We will provide such a comparison in a more 
detailed presentation describing the present work. 

We thank John Baker, Gaurav Khanna, and Raymond 
Puzio for helpful suggestions. This work was partially 
supported by NSF grant PHY-9734871. 
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FIG. 3. Logarithmic representation of the radiation wave- 
forms ReQ (solid curve) and ImQ (dashed curve) as a function 
of t/M at r*/M = 25,6 = tt/4, for a collision with a/M = 0.6, 
and Zo/M = 0.05. See text for discussion. 




FIG. 1. Initial data for ImdtQ as a function of the spatial 
coordinates r* and 6. Note the antisymmetry with respect to 
the equatorial plane, given by 8 — n/2. 




FIG. 2. The relative deviation 5 = ($ — &Keiv)/& is eval- 
uated at f — 0.5, 9 — 0. For a given value of z we compute 
the relative deviation S Z(I , at zq, and the relative deviation 
&2z at 2z . The exponent 7 is defined by 82 ZO = 2 7 c5 Z0 , so 
that 7 is roughly an exponent in the relationship S oc z^ . The 
figure shows that the dependence varies from approximately 
quadratic for small a to approximately linear for large a. 
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FIG. 4. The calculated total energy -E ra d radiated in grav- 
itational waves, as a function of Zo/M for different values of 
a. A meaningful demonstration of the influence of o would 
require a physical measure of separation in place of zq. See 
text for details. 
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